Motion of discrete solitons assisted by nonlinearity management 
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We demonstrate that periodic modulation of the nonlinearity coefficient in the discrete nonlinear 
Schrodinger (DNLS) equation can strongly facilitate creation of traveling solitons in the lattice. We 
predict this possibility in an analytical form, and test it in direct simulations. Systematic simulations 
reveal several generic dynamical regimes, depending on the amplitude and frequency of the time 
modulation, and on initial thrust which sets the soliton in motion. These regimes include irregular 
motion, regular motion of a decaying soliton, and regular motion of a stable one. The motion may 
occur in both the straight and reverse directions, relative to the initial thrust. In the case of stable 
motion, extremely long simulations in a lattice with periodic boundary conditions demonstrate that 
the soliton keeps moving as long as we can monitor without any visible loss. Velocities of moving 
stable solitons are in good agreement with the analytical prediction, which is based on requiring 
a resonance between the ac drive and motion of the soliton through the periodic potential. All 
the generic dynamical regimes are mapped in the model's parameter space. Collisions between 
moving stable solitons are briefly investigated too, with a conclusion that two different outcomes 
are possible: elastic bounce, or bounce with mass transfer from one soliton to the other. The model 
can be realized experimentally in a Bose-Einstein condensate trapped in a deep optical lattice. 



I. INTRODUCTION 



The discrete nonlinear Schrodinger (DNLS) equation is a well-known fundamental model of nonlinear lattice dynam- 
ics, which allows to study many features of the nonintegrable dynamics in the universal setting 0]. Simultaneously, 
this model finds direct applications to arrays of nonlinear optical waveguides (as it was predicted long ago 2] and 
demonstrated in detail more recently, see Refs. and references therein), and to arrays of droplets in Bosc-Einstcin 
condensates (BECs) trapped in a very deep optical lattice 3 . 

In all these contexts, discrete solitons are the most fundamental dynamical excitations supported by the DNLS 
equation. The dynamics of standing solitons, pinned by the underlying lattice, is understood quite well, in terms 
of both numerical simulations and analytical approximations, the most general one being based on the variational 
method 0. However, moving discrete solitons is a more complex issue [1,[3- While, strictly speaking, exact solutions 
for moving solitons should not exist because of the radiation loss, direct simulations indicate that a soliton may move 
freely if its "mass" norm) does not exceed a certain critical value 0. In the quasi-continuum approximation, the 
source of the braking force acting on the moving soliton is the effective Peierls-Nabarro (PN) potential induced by 
the lattice 0- 

In the case of the DNLS equation describing an array of nearly isolated droplets of a BEG in a deep optical lattice, 
an interesting possibility is to apply the ac Feshhach-resonance management (FRM) to it, as it was recently proposed 
in Ref. FRM may be induced by an external ac magnetic field, which periodically (in time) changes the sign of 
the nonlinearity by dint of the FR affecting collisions between atoms (for a one-dimensional BEG without the optical 
lattice, the concept of FRM was proposed in Ref. ) . 

In this work, our aim is to demonstrate that the FRM, applied to the DNLS model, can strongly facilitate the 
motion of discrete solitons (Ref. p| was only dealing with standing solitons). The model is based on the equation 



m„ + W„+i + Un-l - 2Un + g(i)l"nP"« = 0: (1) 

where u„(i) are the mean- field BEG wave functions at the lattice sites, and the real time-dependent nonlinear coeffi- 
cient, proportional to the scattering length, is 

9{t) = 9dc + ffac sin(wi). (2) 

Our presentation will be structured as follows: In section II, approximating the soliton by a Gaussian, we present an 
analytical estimate of the effective PN potential for the moving discrete soliton. The estimate suggests that the ac 
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modulation of g{t) may indeed help to suppress the PN potential, and thus facilitate free motion of discrete solitons. 
In section III, we display results of systematic simulations, summarized in the form of diagrams in the parameter 
plane {uj,gac)- The diagrams feature several generic dynamical regimes, including a large area of stable progressive 
motion, that can last indefinitely long. Collisions between solitons moving in opposite directions are briefly considered 
too, with a conclusion that they bounce from each other, sometimes featuring mass transfer between the solitons. 
Diverse dynamical regimes predicted in this work suggest straightforward possibilities for new experiments in the 
BEC trapped in a deep optical lattice. Finally, in section IV, we summarize our findings and present our conclusions. 



II. ANALYTICAL APPROXIMATION 



The continuum limit suggests the following ansatz for a moving soliton 



Wans(",i) = Aexp 



(3) 



where A,a,^{t), and (j) are, respectively, the amplitude, squared inverse width, central coordinate, and phase of the 
soliton. Accordingly, ^ is the soliton's velocity, simultaneously being the wavenumber of the wave field carrying 
the moving soliton. For the true soliton, in the continuum limit, the variational approximation yields the following 
relations: 



3a, = 4:V2a/g, 



(4) 



if g =const > 0. In Eqs. a is regarded as an arbitrary positive constant (intrinsic parameter of the soliton family). 

To estimate the PN potential acting on the quasi- continuum soliton in the discrete system, we use the Hamiltonian 
corresponding to the DNLS equation (0), 
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(5) 



the asterisk standing for the complex conjugation. Substituting the ansatz Q into H, the potential of the soliton- 

lattice interaction can be calculated in the form of a Fourier scries, H{£^,£^) = X]m=o ^™(^) '^'^^ (^'''"^O- the case 
of a broad soliton, for which the ansatz Q is relevant, we keep only the lowest harmonic (m = 1), which is actually 
the PN potential C/pN- After some simple algebra, we thus find 



c/pN(e,0 
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cos (27r^) , 



(6) 



A novel feature in this estimate, in comparison with known perturbative results is the dependence on the soliton's 
velocity, 

If the soliton's relation between a and A^, as given by Eq. (for constant g), is substituted in Eq. the 
coefficient in front of cos (27r^), i.e., the amplitude of the PN potential, never vanishes. However, it may vanish if the 
underlying pulse ^ is considered not as a soliton, but just as a pulse with A^ and a taken independently; then, the 
condition of the vanishing of the PN potential determines a discrete spectrum of the velocities ^, in the form 



1 + exp (-a/2) cos = (gA^/A) ^ir/ (2a) exp (7rV(4a)) 



(7) 



provided that gA^ is small enough to make the right-hand side of Eq. ^ smaller than 2 [this caveat is essential, as 
the factor exp (7r^/(4a)) may be exponentially large]. 

In this work, however, our objective is not to verify this possibility, but rather to consider the case when the 
nonlinear coefficient g is a function of time, as per Eq. ||2J). Note that, for a broad soliton (small a), the PN potential 
barrier is exponentially small, hence the soliton's kinetic energy may be much larger than the potential. This implies 
that the velocity of the soliton moving through the potential © with the period L = 1 contains a constant (dc) part 
and a small ac correction to it, with the frequency 27r^o/L = 2n^o IM- 



iit) ^io+ii cos (27r4<) , ^2 « 



(8) 
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Then, the substitution of the expression ((SJ into the condition which provides for the suppression of the PN 
potential, one can expand its left-hand side, 

1 + exp (-a/2) cos « 1 + exp (-a/2) cos (^^q /2) - (Ci/2) sin (eo/2) cos (2niot)] • (9) 

Now, substituting the variable g{t) from Eq. (0) into the right-hand side of Eq. ((T)), it is obvious that g^c and gac can 
be chosen so as to provide for the fulfilment of the condition, provided that the average soliton's velocity takes the 
resonant value, = ll!/2tt. More generally, due to anharmonic effects, one may expect the existence of a spectrum of 
resonant velocities, 

^0 = (Crcs)i^''^ = Mlo/27tN, (10) 

with integers M and N. 

Actually, an ac drive can support stable progressive motion of solitons at the resonant velocities l|l()(l (assuming the 
spatial period L = 1), even in the presence of dissipation, in a broad class of systems. This effect was first predicted for 
discrete systems (of the Toda-lattice and Frenkel-Kontorova types) in Refs. [l3l| , and demonstrated experimentally in 
an LC electric transmission line in Ref. ■ Later, the same effect was predicted jTEj and demonstrated experimentally 
|l6l | in continuous long Josephson junctions with a spatially periodic inhomogeneity. All of the above examples were 
in the context topological (kink- like) excitations in the aforementioned models. The only example similar to what is 
suggested here in the context of non-topological solitons that we are aware of, was in the work of 17] (the context in 
the latter case involved damping and external, fixed and localized in space ac-drive). Note also that the mechanism 
of the ac-driven motion considered here is different from that in ratchet systems (see, e.g., Ref. |18j and references 
therein) . 



III. NUMERICAL RESULTS 



We now proceed to examine the above analytical predictions through direct numerical simulations. Wc integrate 
Eq. with an initial configuration in the form of a standing-soliton solution for the case of ^ac = 0. This solution 
has the form of M„(i) = f„ exp(wi), with the real field d„ obeying the equation 

lyVn = Vn+l + Vn-l - 2Vn + gdc^n- (H) 

Equation (|ll|l was solved by means of well-known methods (starting from the anti-continuum limit) . Then, to set the 
soliton in motion, it was provided with a lattice momentum q, so that the initial condition was 

u„(0) = exp(m9/2). (12) 

Equation l|12() implies that the soliton will move to the right if q > 0. 

The results will be displayed for ^dc = 1, = 1, and three values of the initial thrust, q — 0.25, q — 0.5 and q — 1, as 
these cases were found to represent a generic situation in the plane of the ac-drive's parameters (w, (jac)- Simulations 
were run in the interval < t < 100 x (2Tr/uj) or longer, by means of the fourth-order Runge-Kutta algorithm, with 
the time step At = 0.002. In most cases, the lattice with 251 sites was used. Edge absorbers were installed by adding 
the loss term ijUn, with 7 = 1, to Eq. Q at the ten sites adjacent to each edge. Besides that, extremely long 
simulations were performed in a longer lattice with periodic boundary conditions, to verify if the motion could last 
for very long times, and also to examine collisions between the solitons, see below. 

If 5ac — 0, the soliton pushed as per Eq. H12|l with q < 0.7 docs not move. Instead, it remains pinned to the lattice, 
with its center oscillating around an equilibrium position. This observation may be explained by the fact that the 
kinetic energy given to the soliton is smaller than the height of the PN potential barrier. 

Several distinct types of the dynamics were observed at g^c > 0, depending on the driving frequency uj and the 
thrust parameter q. First, the soliton may remain pinned, as shown in Fig. ^ In this case, the simulations [run in the 
interval < t < 300 x (27r/a;)] demonstrate that the soliton stays pinned within a few sites from its initial position. 
The central coordinate and "participation number" P (actually, it measures an average width of the soliton) , shown 
in Fig. Hand elsewhere, are defined as 

x^^nKp/^Kp, p= f^Kp] (13) 
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The next generic regime is that of irregular motion, as shown in Fig. [5] A characteristic feature of this regime 
is that the sohton randomly changes the direction of motion several times, and the velocity remains very small in 
comparison with regimes of "true motion" , see below. 

Under the action of a stronger drive, the soliton can also split into two secondary ones moving in opposite directions, 
see Fig. |21 As is seen, the splitting is strongly asymmetric, and the heavier secondary soliton (splinter) may move 
both forward and backward, with respect to the initial push. We notice that splitting of a quiescent soliton (without 
any initial thrust applied to it) into two symmetric splinters, moving in opposite directions with equal velocities, was 
reported, in the same model based on Eq. (Q, in Ref. ,9]. The initial push applied to the soliton is a natural cause 
for the symmetry breaking observed here in the case of the splitting. 

Two distinct regimes of regular motion of the soliton were also observed. In one case, shown in Figs. 01 and |S1 
the moving soliton is not really stable, as it gradually decays into radiation. The most interesting case is that of 
persistent stable motion of the soliton, without any observable decay (after an initial transient stage of the evolution). 
Examples of the latter are displayed in Figs. |Sland|7| To distinguish between the two different regimes of free motion, 
we have adopted a criterion that the moving soliton is stable if it keeps more than 70% of the mass, X] I"""! i its 
core. While this criterion entails a degree of arbitrariness, we have found that it accurately represents the dynamics 
of the structures resulting from the FRM. The generic examples displayed in Figs. ^-[Tldemonstrate that free motion 
of the soliton is possible in both the straight and reverse directions, relative to the initial thrust. 

As stable motion of solitons in the nonintegrable discrete model is an issue of obvious interest, we have further 
investigated this case, replacing the finite lattice with edge absorbers by a ring-shaped one, with periodic boundary 
conditions. This setting opens a way to study indefinitely long motion of the soliton. The result, illustrated by 
examples shown in Fig. |S1 is that the moving solitons remains stable (preserving its shape) as long as the simulations 
could be run. In this case, it is natural to compare the average velocity c of the persistent motion with the prediction 
given by Eq. H1U|I . The result is ci « 0.246 and C2 ~ 0.155 in the cases shown in Figs. Eta) and and Fig. |Htb), 
respectively. Comparison with the analytical formula (|10|l (with = 1, which is the driving frequency in the examples 
shown in Fig. [HJ demonstrate that C2 and ci fit well to the predicted values in the cases of the, respectively, fundamental 
and second-order resonance, 

C2/ (c.cs)l'^ ~ 0.974, ci/ (c,.es)f ~ 1-029 (14) 

Relatively small discrepancies between the predicted and observed values in Eq. 11411 can be explained by the fact 
that the effective perturbations are not really weak in these cases. 

For the fixed values of gdc = 1 in Eq. ||2Jl and i/ = 1 in Eq. H12I) . and several values of the initial thrust in Eq. 
(initial), q — 0.25, 0.5, and 1, we have collected results of systematic simulations, varying the drive's parameters, (^ac 
and Lu, by small steps in broad ranges spanning the relevant two-parameter space. The results are encompassed in 
Fig. 1^1 in the form of maps in the (w,(7ac) plane, where we outline regions giving rise to each of the qualitatively 
different dynamical regimes described above. This road map is one of the key findings of the present work, detailing 
the various possibilities arising as a result of FRM in a dynamical lattice. 

Some general features can be deduced from examination of the maps in Fig. O As is seen, the increase of the 
initial thrust q significantly affects the map, although quantitatively, rather than qualitatively. At all values of q, 
the irregular motion is, generally, changed by stable progressive motion (straight or reverse) with the increase of 
the drive's amplitude, and/or decrease of its frequency, which seems quite natural. Further increase of the drive's 
strength, which implies the applications of a relatively strong perturbation to the system, may be expected to lead 
to an instability, which indeed happens, in the form of an onset of the gradual decay of the moving solitons. Finally, 
strong instability sets in, manifesting itself in the splitting of the soliton. It also seems natural that the soliton is 
more prone to splitting if the driving frequency is low, as internal strain in the pulse, which eventually leads to its 
splitting, has more time to accumulate if the drive oscillates slowly. 

Transition to the reverse motion tends to happen parallel to the transition from the stable moving soliton to the 
decaying one. For this reason, in most cases (but not always) backward-moving soliton are decaying ones. Finally, it 
should be noted that the increase of the initial thrust leads to overall stabilization of the soliton (somewhat counter- 
intuitively), making the decay and splitting zones smaller. 

Finally, using the large lattice with the periodic boundary conditions, we also simulated collisions between solitons 
originally moving with opposite velocities. Initial pulses were generated by applying the thrust ±q to two quiescent 
solitons. A systematic study of collisions is very difficult in the present model, cf. Ref. |7|]. Nevertheless, we were 
able to identify two different types of the interaction, typical examples of which are shown in in Fig. 1101 

In the case of Fig. Iior a). the solitons bounce back from each other quite elastically. Afterwards, one of the soliton 
spontaneously reverses its direction of motion, due to its interaction with the underlying lattice. Eventually, we 
observe a pair of virtually noninteracting solitons traveling indefinitely long in the same direction. 

In another case. Fig. llOf b). the solitons also bounce after the first collision; however, in this case the collisions are 
inelastic, resulting in transfer of mass from one soliton to the other. Repeated collisions lead to additional transfer, 
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and eventually the weak soliton almost disappears. 



IV. COCLUSIONS 

In this work, we have investigated moving solitons in the DNLS equation with a periodically time-modulated 
nonlinear coefficient. An approximate analytical consideration predicts that the ac nonlinearity management may 
facilitate creation of traveling solitons in the lattice. Systematic simulations reveal several generic dynamical regimes, 
depending on parameters of the time modulation, and the initial thrust which sets the soliton in motion. Besides 
the possibility that the soliton remains pinned, the basic dynamical regimes feature irregular motion, regular motion 
of a decaying soliton, and regular motion of a stable one. In the latter case, extremely long simulations in the 
lattice with periodic boundary conditions demonstrate that the soliton keeps moving indefinitely long without any 
tangible loss. Velocities of the moving stable solitons are found to be in good agreement with the analytical prediction 
through a resonance condition, revealing the key mechanism for sustaining the travelling motion. All the generic 
dynamical regimes were mapped in the model's parameter space. Collisions between stable moving solitons were 
briefiy investigated too, with a conclusion that two different outcomes are possible: elastic bounce, and bounce with 
mass transfer from one soliton to the other. The model can be realized in a Bose-Einstein condensate trapped in a 
deep optical lattice, with the nonlinearity modulation induced by the Feshbach resonance in the ac regime. 

It would be interesting to examine similar dynamical scenaria for the effect of FRM on dark solitons, and also to 
monitor how this phenomenology is modified in higher dimensions. Such studies are currently in progress and will be 
reported in future publications. 
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(a) (b) (c) 




FIG. 1: A typical example of a dynamical regime in which the soliton remains pinned, for g^c = 0.03, = 1 and q — 0.5. Here 
and in similar plots below, panels display the following variables as functions of time: (a) positions of the center of mass X 
(red) and maximum Xp (blue) of the density field, |?inp; (b) the "participation number" P (which characterizes the width of 
the soliton); (c) the density amplitude of the soliton, i.e., \u„f at n = Xp. 

(a) (b) (c) 



'IP 

'(lllllli!,!!. I , - 

100 200 300 400 500 600 100 200 300 400 500 600 100 200 300 400 500 600 

t t 1 




FIG. 2: A generic example of irregular motion of the soliton, for ^ac = 0.065, u = 1 and q = 0.5. 
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FIG. 3: An example of asymmetric splitting of the soliton, for ~ 0.196, ui = 0.5 and q = 0.5. Unlike the previous figures, here 
panels (a) and (c) show, respectively, the position and magnitude of two local density maxima corresponding to the secondary 
solitons (splinters) past the splitting point, and panel (b) shows, accordingly, two "participation numbers", computed as per 
Eq. 1131 . but separately for n < and n > 0. In this figure and below, an additional panel (d) shows the initial and final 
density distributions, \u„{t — 0)f and \u„{t = ten|^, respectively, by the dashed and continuous lines. In this case, tun — 300. 
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(a) (b) 




FIG. 4; A generic example of the free motion of a decaying soliton in tiie straight direction, for gac ~ 0.206, uj — 0.5 and 
q = 0.5. In panel (d), the final configuration pertains to tun = 400. 




FIG. 5: The same as in Fig. 0] in the case of motion in the reverse direction, for gac ~ 0.218, lu = 0.5 and q — 0.5. In panel 
(d), the final configuration is shown for tun = 350. 
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FIG. 6: A generic example of free motion of a stable (nondecaying) soliton in the straight direction, for gac ~ 0.132, ui — 1 and 
q = 0.5. The panel (d) shows the final configuration at tun = 500. 




FIG. 7: The same as in Fig. |S| but in the case of motion of a stable soliton in the reverse direction, for pac ~ 0.122, cu = 1 and 
q = 0.5. The panel (d) shows the final configuration at tfln = 500. 




FIG. 8: Indefinitely long motion of stable solitons in a ring lattice composed of 601 sites. The gray-scale plots show the 
spatiotemporal distribution of the density, |Mn (t) l'^ . The examples of the straight (a) and inverse (b) motion pertain, respectively, 
to 3ac = 0.132, uj^l,q^ 0.5, and g^c = 0.122, cu ^l,q = 0.5. 
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q = 0.25 
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FIG. 9: Maps in the left column show areas in the parameter plane {uj,gax:) which give rise to the following dynamical regimes. 
White areas: the soliton remains pinned; cyan: irregular motion; green: splitting; magenta: regular motion with decay; black: 
stable motion (without decay). The maps in the right column additionally show the difference between the forward (alias 
straight, marked by magenta) and backward (alias reverse, marked by black) directions of the regular motion, relative to the 
direction singled out by the initial push. Regular-motion regimes for both decaying and stable solitons are included here. 
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FIG. 10: Two typical examples of different outcomes of collisions between solitons with equal masses moving in opposite 
directions in the lattice with periodic boundary conditions. The parameters are (jac = 0.132, ui = 1, q = 0.5 (a) and gac = 0.122, 
uj = 1, q = 0.5 (b). Note that the collision is multiple in panel (b). 



